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Abstract: The DAMA collaboration have claimed to detect particle dark matter (DM) 
via an annual modulation in their observed recoil event rate. This appears to be in strong 
disagreement with the null results of other experiments if interpreted in terms of elastic 
DM scattering, while agreement for a small region of parameter space is possible for in- 
elastic DM (iDM) due to the altered kinematics of the collision. To date most analyses 
assume a simple galactic halo DM velocity distribution, the Standard Halo Model, but di- 
rect experimental support for the SHM is severely lacking and theoretical studies indicate 
possible significant differences. We investigate the dependence of DAMA and the other 
direct detection experiments on the local DM velocity distribution, utilizing the results of 
the Via Lactea and Dark Disc numerical simulations. We also investigate effects of varying 
the solar circular velocity, the DM escape velocity, and the DAMA quenching factor within 
experimental limits. Our data set includes the latest ZEPLIN-III results, as well as full 
publicly available data sets. Due to the more sensitive dependence of the inelastic cross 
section on the velocity distribution, we find that with Via Lactea the DAMA results can be 
consistent with all other experiments over an enlarged region of iDM parameter space, with 
higher mass particles being preferred, while Dark Disc does not lead to an improvement. 
A definitive test of DAMA for iDM requires heavy element detectors. 
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1. Introduction 



First-hand information on the nature of the dark matter (DM) in our universe may be 
revealed by one or more of the current generation of direct detection experiments. The 
DAMA collaboration have measured an annual modulation in their scintillation event rate, 
initially with the DAMA/Nal set-up [1], and more recently with the DAMA/LIBRA con- 
figuration [2], which they interpret as being due to the modulation in the DM-nuclear 
scattering event rate following from the annual variation of the Earth's velocity through 
the DM in our galaxy. At face value, this measurement appears to be in disagreement with 
other direct detection experiments if we make the usual assumptions about the nature 
of the DM; namely that the DM is a weakly interacting massive particle (WIMP) that 
elastically scatters off the target material with a spin independent interaction. 

Inelastic dark matter (iDM) [3] was proposed as a way to reconcile the positive re- 
sult from DAMA/Nal with the null result from germanium based detectors. In the iDM 
scenario, WIMP-nucleon elastic scattering is suppressed, while inelastic scattering from a 
ground-state WIMP to a slightly higher mass excited WIMP is allowed and dominates the 
recoil event rate. As we summarize in Section 2, the kinematics of the recoil scattering 
are changed by the inelastic nature of the collision, and this can bring the DAMA results 
closer to agreement with the other experiments for three principal reasons: 

• Heavier nuclei are favoured. 

• The recoil spectrum is changed at low energies. 

• The ratio of modulated to unmodulated signal is higher. 

Nevertheless the other direct detection experiments still strongly constrain the DM inter- 
pretation of DAMA/LIBRA, with, apparently, only a relatively small region of parameter 
space being allowed [4]. 

One unchecked assumption that goes into the analysis is that the velocity distribution 
of the DM in the galactic rest frame is well described by so-called Standard Halo Model 
(SHM) . The SHM assumes that in the galactic frame the WIMP distribution is an isotropic 
isothermal sphere, which leads to an essentially structureless isotropic Maxwell-Boltzmann 
velocity distribution with dispersion set by the local circular velocity. This is a questionable 
assumption: Little is known about the phase-space distribution of DM on the scales relevant 
for direct detection experiments, and as we outline in Section 3 numerical simulations of 
DM distributions for Milky- Way-like galaxies lead to results differing from the SHM in 
potentially significant ways, especially for iDM which is more sensitive to the DM velocity 
distribution. 

In this paper, we compare the results obtained using the SHM to those obtained using 
two recent computer simulations of the DM distribution in a Milky- Way-like galaxy: 

• Via Lactea [5] - a Milky Way size DM halo distribution containing 234 million par- 
ticles. 
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• Dark Disc [6] - a simulation which contains in addition to the SHM, a slowly rotating 
disc of DM 

With these less idealised halo models, we find that the iDM scenario allows the DAMA 
region to be consistent with all other experiments, and that in the Via Lactea distribution, 
more parameter space can be opened up at high WIMP masses, compared to the SHM. 
After reviewing the details of each experiment and describing how we calculate the allowed 
DAMA region and exclusion curves for the null experiments in Section 4, we present our 
results together with their physical interpretation in Section 5. 



2. Inelastic dark matter 



2.1 Review of inelastic dark matter 



If the dark matter particle can only scatter off nuclei by making a transition to a heav- 
ier state then the altered kinematics of the interaction can lead to significant changes in 
detection rates for different detectors [3, 4, 7]. If we call the two dark matter states X- 
and x+1 with mass splitting 5 = M^-^ — ~ 0(100 keV), then the minimum velocity 
to scatter off a nucleus and impart an energy is given by; 



1 (M^En^^. ^2.1) 



2MnEr V A* 

where // is the reduced mass of the WIMP-nucleus system, and Mjy is the nucleus mass. 
There are a variety of particle physics models that lead to such phenomenology [3, 8, 9] . The 
models will not be our concern in this paper, rather our focus will be upon the consequences 
of iDM for direct detection experiments. 

The differential event rate for WIMP nucleus scattering on a given element, as a 
function of recoil energy Eji, is given by, [10]: 

- iV,M„^(M±A(A^f f (2.2) 



dEn 2m^nl^ /„2 Jv^,4En 

where, as in [4], Nt is the number of target nuclei, is the local dark matter density, 
taken as 0.3 GeV/cm^, fine the reduced mass of the WIMP and nucleon, and (T.„ is the 
WIMP-neutron cross section. The factors /„ and fp parameterize the relative scattering 
strength off neutrons and protons. For simplicity and to maintain model-independence we 
assume that scattering off neutrons and protons is the same and take fn = fp = 1- (This 
assumption is not always correct for candidate dark matter particles from specific models.) 
F[Epi] is the nuclear form factor, which can be modeled analytically with the Helm form 
factor [11] or calculated numerically using the Fermi Two-Parameter form factor [12]. In 
accordance with [4] we wih use the Helm form factor for CDMS II, KIMS, XENONIO, 
ZEPLIN-II, and ZEPLIN-III, and the Fermi Two-Parameter form factor for CRESST-II 
and DAMA/LIBRA. In the first panel of Figure 1 we show for tungsten the difference that 
results between the form factors. 
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Finally, f{v), is the local velocity distribution of the WIMPs in the Earth rest frame. 
The default assumption for the dark matter phase space density is the so-called "Standard 
Halo Model" (SHM) which assumes a local isothermal and isotropic distribution of dark 
matter leading to a Maxwell-Boltzmann velocity distribution. Changes in this velocity 
distribution make major differences to the final differential event rate as calculated from 
Eq.(2.2), especially for iDM. We emphasize that f{v) is poorly constrained by data and not 
well understood. The primary subject of this paper is to investigate the allowed parameter 
space for iDM taking account of reasonable variations in f(v) motivated by numerical 
simulations of the DM distribution in our galaxy. 

Tungsten Spectrum, 6 = 100 keV Tungsten Spectrum, 6 = keV 

Rate (cpd/ka/keV) Rate (cpd/kg/keV) 

0.04r " 5r 




Figure 1: Recoil energy spectrum for scattering off of tungsten. On the left the inelastic scattering 
rates are shown for the two choices of form factor: the Fermi Two-Parameter (dashed) and Helm 
(solid) form factor. Both choices show that iDM leads to a suppression of low-energy events. The 
right panel shows the recoil energy spectrum for clastic scattering; in this case the difference between 
the two form factors is negligible. For elastic scattering the recoil spectrum peaks at low En. All 
calculations assume M-^ = 200 GeV, (t„ — 10^^° cm^ and the SHM velocity distribution truncated 
at Wesc — 500 km/s. 

The most important difference between inelastic and standard elastic dark matter is 
the (5-dependent increase in the minimum relative velocity for scattering with recoil energy 
Er, see Eq.(2.1). The important consequences for direct detection experiments are: 

• As only higher velocity portions of the dark matter distribution lead to recoils the 
overall event rate for a given cross section is lower. 

• The spectrum of events is changed, with the greatest qualitative departure from 
elastic scattering occurring at low energies. As one can see from the first panel of 
Figure 1 there is a low energy cut off in the recoil spectrum. This is not present for 
elastic scattering as shown in the second panel of Figure 1. 

• Due to the higher minimum velocity the annual modulation as a fraction of the 
total signal can be much greater. It is possible that detection of WIMPs occurs 
solely during the summer months, when the Earth's velocity into the WIMP wind is 
greatest, leading to a modulation fraction of 100%. This has important consequences 
both for DAMA/LIBRA and XENONIO. 
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• As the minimum velocity cut-off, Eq.(2.1), is lower for heavier target nuclei, inelas- 
tic scattering leads to higher expected detection rates for heavier elements. This 
has important consequences for the sensitivity of the CDMS II {Ace = 72.64) and 
CRESST-II {Aw = 183.84) experiments. 

2.2 Overview of consequences of iDM for the various experiments 

The minimum velocity for scattering will be much greater for the CDMS II detector com- 
pared to DAMA/LIBRA, where the dominant signal is from scattering off of germanium 
and iodine respectively, because germanium is much lighter than iodine. This leads to a 
weakening of the constraints from CDMS II on the region of parameter space preferred by 
the DAMA/LIBRA annual modulation signal. In fact for high enough 5 the expected rate 
for CDMS II can be consistent with zero, while still allowing signal at DAMA/LIBRA. 
The converse effect applies to the CRESST-II experiment, where scattering occurs off of 
tungsten. Here we would expect a higher rate for scattering, which leads to CRESST-II 
setting the most stringent constraints on the DAMA/LIBRA results. 

The three xenon based experiments, XENONIO, ZEPLIN-II, and its successor ZEPLIN- 
III, provide an excellent model-independent test of the DAMA/LIBRA results due to the 
similarity of xenon and iodine masses {Axe = 131.293,^/ = 126.904). However the events 
observed at these detectors lead to consistency with the DAMA/LIBRA preferred region. 
Also, due to the enhancement of the modulation effect, detectors which took results over 
the winter period, when signals would be lowest, inevitably set lower constraints than pos- 
sible if running during the summer months. This applies in particular to the XENONIO 
experiment. Moreover, the low upper-energy limit (30.2 keV) in the analyzed recoil energy 
spectrum of ZEPLIN-III reduces the sensitivity of this experiment to the iDM scenario, as 
for a typical 5 of 100 keV the xenon recoil spectrum peaks at ~40 keV. 

We return to the details of the experiments in section 4. 

3. Dark matter halos 

Due to the long range of the gravitational force, and the distribution of matter in our galaxy, 
one would expect the correct velocity distribution of dark matter particles to deviate from 
exact Maxwellian and to show some anisotropy. N-body simulations of large numbers of 
dark matter particles have shown that the SHM may well be incorrect [5, 6]. 

Here we will consider three models for the local velocity distribution of particles in 
the dark matter halo. We will always take the Earth velocity with respect to the galactic 
rest frame to be given by v Earth = 1^0 + '^^ei where Vq is the Sun's velocity relative to the 
galactic rest frame and is the Earth's velocity relative to the Sun. Vq is the sum of the 
Sun's peculiar [13] and circular [14] velocities: 




(3.1) 



V(Q is given by [10] ; 
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■u© = {ue){1 - esin(A(n) - Aq)) 



^ cos{l3x) sin(A(n) - 

cos(/3y) sin(A(n) — Aj,) | km/s (3.2) 
^ cos{l3z) sin(A(n) - 



where the Earth's orbit has a mean velocity {ue) = 29.79 km/s and eUipticity e = 0.016722. 
The quantities Pi, Aj, define the orientation of the Earth's orbit in galactic coordinates, 
and A(n) gives the angular position of the Earth's orbit for a given day number n, with 
n = 1 corresponding to 1st January 2000. These quantities are given in [10]. 

The WIMP velocity distribution in the Earth's rest frame is found by performing a 
Galilean boost from the galactic rest frame distribution: 

/e(f , t) = fgal{v + VEarth{t)) (3.3) 

One feature common to all of the velocity distributions considered here is that they are 
truncated at the local escape velocity, fgai{\v\ > Vesc) = 0. There is relatively large error 
in the known value of the local escape velocity, 498 < v^sc < 608 km/s (90%), with a 
median of 544 km/s [15]. We have taken Vesc = 550 km/s in the analysis presented here. In 
Section 5.6 we check that varying Vesc between the 90% confidence limits does not change 
the qualitative features of our results, and therefore the conclusions. 

3.1 The Standard Halo Model 

In the Standard Halo Model (SHM) one assumes a Maxwellian velocity distribution for 

the dark matter particles, with a local velocity dispersion a = v/\/^, where v = \J 

is the local circular velocity, and U (r) is the gravitational potential. This distribution is 

truncated at the escape velocity according to: 



V > Vp 



3.2 Via Lactea 



Here we use results for the phase space distribution of dark matter in a Milky- Way-like 
galaxy derived from a simulation containing 234 million particles of dark matter and no 
baryons; Via Lactea [5] published in 2006. 

In [16] the velocity distribution of dark matter particles was fitted to the distribution of 
individual particles from the Via Lactea simulation [5]. The radial and tangential velocity 
distributions were fitted according to; 



ff \ ^^^^ 



4 



Or- 

(3.5) 



aj'- 



(3.6) 



and the distribution was truncated in the same way as for the SHM. From Figure 3 in 
[16] we extracted the values ^ 1.09, and vr/{\/—U {vq)) « 0.72 and ay ~ 0.73, and 
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vt / {\/ —Uiro)) 0.47 at our radius from the centre of the galaxy, tq f« 8.5 kpc. (Although 
there is uncertainty in our radius from the galactic centre, tq = 8.0 zt 0.5 kpc [14], from 
Figure 3 in [16] one can see that the values of ut and vt do not change significantly over the 
range 7.5 < tq < 8.5 kpc, and the tangential velocity distribution is the main determinant 
of the event rate) . 

For completeness we will present results using two values for f/(ro). In the Via 
Lactea halo the average value of \/—U {r) between 7 and 9 kpc is 270km/s [17] ; we will refer 
to results using this value as VL270. Following [16], and in order to allow direct comparison 
between inelastic and elastic scenarios, we also present results using y/-U{ro) = 220km/s 
[18], and refer to these as VL220- 

We believe both values for \/—lJ{ro) are worth studying, given that the Milky Way is 
baryon dominated at the solar radius, so uncertainties will arise from the lack of baryons in 
any simulation which contains only dark matter particles. It is also important to remember 
that the velocity distribution extracted from any DM simulation is for a Milky- Way-like 
galaxy, and not the Milky- Way itself. Therefore it is important to study the effects of 
reasonable deviations from results predicted by a simulation. 

The deviation from the Gaussian distribution and difference in radial and tangential 
velocity dispersions have been shown in [16] to affect the expected DAMA/LIBRA mod- 
ulation signal for elastic dark matter, although not by enough to allow an elastic DM 
interpretation of the DAMA results. 

Due to the smaller tangential velocity dispersion in this model, compared to the SHM, 
one would expect the results for inelastic scattering to be changed for light nuclei. This 
is because the high minimum velocity for scattering on light nuclei leaves only the high 
velocity part of the distribution detectable. This smaller velocity dispersion reduces the 
population of this high velocity region further, leading to a reduced event rate. This effect 
is particularly interesting as it arises through the combination of iDM and the Via Lactea 
halo. We further discuss this in Section 5.1. 

3.3 Dark Disc 

Previous simulations of the DM phase space density distribution have modelled the dark 
matter alone, while at the solar neighbourhood we expect the effects of the baryons, the gas 
and stars that make up the Milky Way, to be important. Read et al. [6] have performed a 
series of simulations including the baryons, and have shown that massive satellites, dragged 
into the disc plane by dynamical friction, are torn apart by tidal forces depositing their 
stars and dark matter into a thick disc, lying in the same plane as the visible galaxy. Later 
the influence of this dark disc on direct detection of elastic dark matter was examined [19]. 
Here we investigate the influence of the dark disc on iDM detection. 

Following [19] , we assume the dark disc kinematics match the Milky Way's stellar thick 
disc, whose properties are listed in Table 1 of [6]. This will be a good approximation if the 
Milky Way's stellar thick disc is mostly composed of accreted, rather than heated stars. 
We model the dark disc as a component of DM additional to the SHM, lagging the rotation 
of the Sun by 40 km/s in the tangential direction, compared to the SHM which lags by 
220 km/s. It is also assumed to have a Maxwellian velocity distribution, with a dispersion 
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of o" = (63,39,39) km/s, and a density in the range 0.5 < < 2. We employ these 

parameter ranges in our study of iDM signals from the dark disc. 

Due to the smaller relative velocity of the dark disc with the Earth, when compared 
with the SHM, one would expect dark disc event rates to be low for high 6 as the dark disc 
particles may not have high enough velocity to cross the minimum velocity threshold. 



4. Experimental information 

We now turn to a detailed discussion of the individual direct detection experiments. 
4.1 DAMA/LIBRA 

As the Earth revolves around the Sun, there should be a larger flux of WIMPs incident on 
the detector around the 2nd June, when the relative velocity of the Earth is at a maximum 
with respect to the galaxy. Conversely, the flux incident on the detector should be smallest 
around the 2nd December when the relative velocity of the Earth is at a minimum with 
respect to the galaxy. It is the annual modulation in the recoil event rate caused by this 
velocity modulation that the DAMA collaboration claim to have measured. 

The first results from the DAMA/LIBRA set-up, with an exposure of 0.53 ton-yr, 
have recently been published [2] . These have been combined with the data collected by the 
DAMA/Nal set-up to give an impressive 0.82 ton-yr total exposure, yielding a modulation 
signal at 8.2 a C.L. 

For a WIMP mass > 10 GeV, this signal appears to be in conflict with other 
experiments under the assumptions of spin independent, elastic WIMP-nucleon scattering, 
with the WIMP phase space distribution described by the SHM. This has prompted a 
number of alternative explanations: light WIMPs [20, 21, 22, 23, 24, 25, 26, 27], spin 
dependent interactions [28], mirror dark matter [29, 30] and iDM [3]. 

4.1.1 Quenching and channeling 

DAMA use highly radiopure Nal(Tl) scintillators as their target material. The light yield 
of scintillators depends on whether the recoiling nucleus interacts electromagnetically or 
via the strong nuclear interaction, since only electromagnetic interactions will produce 
photons. As a result, the measured energy, Em, is different from the recoiling energy of 
the nucleus. This difference is expressed by the quenching factor q defined by Em = qEr. 
We will follow the convention of measuring Er in keV, and Em in keVee (keV electron 
equivalent). The DAMA collaboration have measured the quenching factor for iodine and 
sodium for their Nal(Tl) crystals. They obtained the values q^a = 0.3 and qi = 0.09 
[31] with an error of 0.01 on the value of qj [32]. Following the analysis of Ref.[4], in our 
calculations, we use qNa = 0.3 and qj = 0.085, however in Section 5.4 we investigate what 
effect varying qj within the experimental limits has on our results. 

If the nucleus recoils along certain directions in a crystalline structure, and if the 
recoiling energy is low, no nuclear interactions occur, so the quenching factor, q, equals 
one. This effect is known as channeling [33]. Since Nal(Tl) is a crystalline material, we 
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need to include this effect in our calculations. We use the parameterisation given in [16] 
for the fraction / of channeled events relevant for DAMA 

^^'^(^^^ " l + 0.75En' ^'^^^^ " 1 + 0.65ER- ^^'^^ 
4.1.2 Calculating limits 

DAMA have released their binned data from the combined DAM A/LIBRA and DAMA/Nal 
data sets covering the range 2-20 keVee. A clear modulation signal is present below about 
8 keVee, while the modulation is consistent with zero at higher energies [2]. The rate in the 
lowest energy bin, covering the range 2-2.5 keVee, is smaller than in the next bin, covering 
the range 2.5-3 keVee, suggesting that the rate is falling to zero at low energies. While 
care should always be taken when examining the data at the edges of the experimental 
sensitivity, it should be noted that both of these features are found naturally with inelastic 
scattering. 

We use a goodness of fit test to analyse the DAMA results. We construct a 
function using the twelve 0.5 keVee width bins between 2.0 - 8.0 keVee and their relative 
uncertainties. We do not fit to the higher energy bins because, as mentioned above, the 
inelastic spectrum falls off at high energies. 

In our 2D plots, we find the best fit point by minimising the function for the 
two unconstrained parameters, either cj„ and M^, or cj„, and 6. For a goodness of fit, 
we require that the best fit point have Xmin < -'-O; given that we have 12 bins and 2 
free parameters. Allowed regions at a given confidence limit are obtained by looking for 
contours = Xmin + ^X^; where A^^ = 4.61 or 10.60 for 90% and 99.5% confidence limits 
respectively. 

DAMA have also released the unmodulated rate for a single-hit scintillation as mea- 
sured by DAMA/LIBRA. A WIMP is expected to scatter once, so we are using the rest 
of the detector as a veto. No other background subtraction is applied. We can use this 
data to set a limit by requiring that the calculated unmodulated rate not be larger than 
the measured rate across the energy range shown in Figure 1 of [2]. 

It should be noted that we have not included the effect of the finite energy resolu- 
tion of the detector in our calculations. We are fitting to the modulated rate from the 
combined data sets of DAMA/LIBRA and DAMA/Nal, however the energy resolution for 
DAMA/LIBRA and DAMA/Nal are known to be different [34, 35]. 

4.2 Null experiments 

In this section, we discuss the data used in calculating the regions excluded at 90% confi- 
dence level by other direct detection experiments. We adopt a conservative approach and 
include events that the experiments ascribe to background processes when setting limits. 
We will find that the most constraining experiment is CRESST-II. Unless otherwise stated, 
we use the "pmax method", described in [36], in setting our limits. This method has the 
advantage of setting similarly strong limits as the "optimum interval method" , but is easier 
to implement. Care should be taken when comparing exclusion curves at 90% confidence 
level and the region allowed by DAMA at 90% confidence limit, since we are using different 
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statistical methods in evaluating them. In calculating the total event rate, we integrate 
the differential rate over the energy range of the experiment, then find the average value 
over the running dates given by the collaboration, and finally, we multiply by the exposure. 
As a final check of our methods, we have reproduced the published exclusion curves from 
CDMS II and XENONIO in the 5 = hmit. 

4.2.1 CDMS II 

In setting the limits for CDMS II, we use the three runs from the Soudan Underground 
Laboratory which were sensitive to nuclear recoil energies between 10-100 keV [37, 38, 39]. 
We only consider scattering from germanium as scattering off of silicon is highly suppressed. 
The published effective germanium exposures are weighted for a WIMP mass of 60 GeV 
and averaged over recoil energies 10-100 keV, however these numbers are expected to apply 
for an inelastic WIMP since their acceptance efficiency is fairly constant over their energy 
range. The first run took place from 11th October 2003 to 11th January 2004, had an 
exposure of 19.4 kg-day and saw one event at 64 keV [37]. The second run was from 25th 
March 2004 to 8th August 2004, had an exposure of 34 kg-day and saw one event at 10.5 
keV [38]. The latest five-tower run was between October 2006 and July 2007, had an 
exposure of 121.3 kg-day and saw no events [39]. 

4.2.2 CRESST-II 

For the CRESST-II limits, we use data collected by the Julia and Daisy detectors [40] 
between 31st January 2004 - 23rd March 2004, and the Verena and Zora detectors [41] 
between 27th March 2007 - 23rd July 2007. We only consider scattering off of the tungsten 
atoms in the CaW04 crystals, and include events on or beneath the curve where 90% of the 
tungsten recoils are expected. The Julia detector had a tungsten exposure of 6.26 kg-day 
and we include the four observed events between 10 - 50 keV and an unpublished event 
above 50 keV [42]. The Daisy detector had a tungsten exposure of 6.84 kg-day and we 
include the two events between 10-12 keV, the event which lies on the curve where 90% of 
the tungsten recoils are expected, at approximately 22 keV, the event at approximately 45 
keV, and a second unpublished event above 50 keV [42]. The Verena and Zora detectors 
had a combined tungsten exposure of 30.6 kg-day and we use the seven published events 
between 10 - 100 keV. To take into account the difference in energy regions, we do not 
integrate above the energy of the extra events for the Daisy and Julia detectors. 

Note that following [4] we have set limits by including the Daisy point which lies on 
the tungsten curve at 22 keV in Figure 9 of Ref.[40]. If we redo the analysis by ignoring 
this point and only including data which lie fully below the curve, we find that CRESST-II 
just manages to exclude all of the DAMA/LIBRA 99.5% allowed region. This point has 
such a large impact on the limits set using the "pmax method" because it lies at the peak 
of the recoil energy spectrum, as can be seen in the left panel of Figure 1. 

This shows us that with further data, the CRESST-II experiment has the ability to 
completely exclude the DAMA/LIBRA region. 



- 10 - 



4.2.3 KIMS 

From the null experiments we consider, KIMS [43] is the only one where scattering occurs 
off of iodine, in the CsI(Tl) crystals they employ. Two of the crystals collected experimental 
data between June 2005 and March 2006, then another two crystals were installed, running 
between December 2005 and March 2006 [44]. In setting limits, we find the average rate 
and statistical error from the four crystals, and require the calculated rate to be less than 
the measured rate plus 1.64 times the error in the first five bins, corresponding to the 
energy range 3-8 keVee. For the quenching factor, we use the parameterization 

_ , 0.175e-^«/i37 
qcsm){En) ^ 1^0 00091^^- (4-2) 

This function fits the curve shown in Figure 13 of [45] to a good accuracy. 

4.2.4 XENONIO 

XENONIO [46] is a liquid xenon based experiment which ran from 6th October 2006 until 
14th February 2007 and had an effective exposure of 316.4 kg-day. In calculating our limits, 
we include the 10 events used in their analysis in the range 4.5 - 26.9 keV, as well as the 14 
events in the range 26.9 - 45 keV. Below 26.9 keV, we use the published acceptances and 
above 26.9 keV we use the known software cut efficiencies and assume a constant nuclear 
recoil acceptance of 0.45 [47]. 

4.2.5 ZEPLIN-II 

We include only the limits from ZEPLIN-II [48] since ZEPLIN-I [49] is not competitive with 
the other experiments we consider [4] . ZEPLIN-II is another xenon based experiment with 
an effective exposure of 225 kg-day and was assumed to run between May and July 2006, as 
suggested by [50] . 29 events were observed in the range 5-20 keVee, which corresponds to 
a range of 13.9 - 55.6 keV if the published quenching factor q = 0.36 is used. In calculating 
our limits, we use the published efficiencies and we take all observed events as signal. 

4.2.6 ZEPLIN-III 

We also include the first results from the ZEPLIN-III experiment [51], the successor to 
ZEPLIN-II. ZEPLIN-III is also xenon based and improves on the limits set by ZEPLIN-II 
as only 7 events were observed, in the range 10.7 - 30.2 keV. The effective exposure was 

126.7 kg-day, the data being collected between 27th February 2008 and 20th May 2008. 
Using Figure 15 of [51], we extract an approximate quenching factor parameterization: 

qxe{EM) ~ (0.142SM + 0.005) exp[-0.305(SM)°-564]_ (4 3) 

We extracted the energies of the 7 observed events from Figure 12 of [51] and converted to 
nuclear recoil energies using the above quenching factor. 

5. Results and discussion 

In this section we present our results for iDM in the context of the different halo models. 
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5.1 VL220 vs SHM 



In Figures 2 and 3 we present a comparison of the iDM limits on the DAMA/LIBRA 
preferred region under the assumptions of the SHM and the VL220 halo, for a variety of 
WIMP masses. 

As one can see the allowed cross section limits are generally higher for the VL220 
halo than the SHM. Also the CRESST-II results are slightly more constraining on the 
DAMA/LIBRA preferred region at low mass when the VL220 halo is used. The most 
interesting feature to note is that all of the experiments show less disagreement with the 
DAMA/LIBRA region at high WIMP masses when the VL220 halo is used. This feature 
is most prominent for the CDMS II exclusion line, but is true for all experiments to some 
degree. 

This feature is made more explicit by observing how the limits change as a function of 
for fixed 5. In Figure 4 we show the exclusion limits for 5 = 100 keV and 50 GeV < 
< 1000 GeV. 

In this case all experiments are less constraining on DAMA/LIBRA at high mass, 
however the CDMS II limits completely rule out the entire DAMA/LIBRA preferred region 
under the assumptions of the SHM, yet for the VL220 halo there is agreement up to ~ 
0( TeV). This can be explained by a combined effect of iDM and the VL220 velocity 
distribution, as noted in Section 3. 

Looking quantitatively at this effect we can take the example of germanium in the 
CDMS II detector. As the inelasticity has pushed the minimum velocity up significantly 
we need only consider the tangential velocity distributions (in the Earth's frame), as they 
are centered around Vdrc ~ 220 km/s whereas the radial velocity distribution is centered 
around vor ~ km/s. Therefore the high velocity components scattering on the germanium 
will be from the tangential component of the halos. 

The CDMS II detector is sensitive to the energy range 10 keV < Er < 100 keV, and 
if we consider 50 GeV < < 500 GeV, then for 6 = 100 keV the minimum velocity 
lies in the range 550 km/s < Vmin < 1006 km/s. As the escape velocity is taken as 
Vesc = 550km/s in the halo rest frame, then the highest velocity particles observed in 
the Earth's frame have Vmax ~ 770km/s. Therefore the observable particles will have 
velocities 550 km/s < Vobs < 770 km/s. In Figure 5 we plot the two different tangential 
velocity distributions and the region over which Vmin varies as is increased. 

As one can see the detectable particle number in the SHM distribution grows signifi- 
cantly as is increased and Vmin decreases, leading to stronger limits on the cross section 
from the CDMS II results. However, for the VL220 distribution this growth in the number 
of detectable particles is not nearly as pronounced, due to the reduced width of the dis- 
tribution, and this leads to much slower variation in CDMS II exclusion limits as is 
varied, and overall to much weaker limits than from the SHM. As stated above this effect 
is due to the combination of iDM, which increases the minimum velocity relative to elastic 
scattering, and the VL220 velocity distribution, which is much smaller than the SHM at 
high velocities. 
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SHM, M. = 80 GeV SHM, My = 100 GeV 
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DAMA Unmodulated 




DAMA Modulated 90% 
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Figure 2: Here we show the variation in the exclusion Hmits set by the experiments as 6 is varied 
and is held constant. These limits are calculated using the SHM. The preferred region of 
parameter space for the DAMA results is shown at 90% and 99.5%, and the DAMA best fit point is 
plotted with a dot. As one can see there is a small region of agreement between all experiments for 
low masses and S - 130 keV. At higher masses both CRESST-H and CDMS H exclude the DAMA 
results and the region of agreement with the other experiments is greatly reduced. 



- 13 - 



VL220, = 80 GeV VL220, = 100 GeV 
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Figure 3: Change in limits when the SHM is replaced by the VL220 halo, cf, Figure 2. At low 
masses there is a smaller region of agreement between CRESST-II and DAMA, with CRESST-II 
almost excluding DAMA at the 90% level over all masses. At high masses, however, CRESST-II 
and CDMS II are significantly less constraining on the DAMA region than for the SHM, this can 
be seen by comparing the bottom right panels of both figures. One can also see that the typical 
cross sections are an order of magnitude higher for the VL220 halo at low masses than for the SHM 
(note also that the lower panels have a different scale for the cross section). 



5.2 VL270 vs SHM 

In Figures 6 and 7, we present the exclusion limits for VL27o- Comparing Figure 6 with 
Figures 2 and 3, and Figure 7 with Figure 4, we see that although the exclusion limits from 
VL270 typically occur at a larger cross section than the SHM, the DAMA allowed region 
is much closer to that from the SHM than VL220- 

With reference to Figure 5, we see that the departure from the SHM velocity distribu- 
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SHM, 6 = 100 keV VL220, S = 100 keV 




[GeV] [GeV] 

Figure 4: The allowed parameter space for fixed 6 = lOOkeV and varying M^. For the SHM (left 
panel) and this value of 6, CDMS II excludes the DAMA region at 90%. For the VL220 halo (right 
panel) the tightest constraints are set by CRESST-II, and there is agreement between DAMA and 
CDMS II up to high WIMP masses. Again one can see that the typical allowed cross sections are 
an order of magnitude higher for the VL220 halo. 

f(V7-) 




Figure 5: Detectable particle distributions in a germanium detector as a function of tangential 
velocity. The detectable region is shaded in grey. The solid, dashed and dot-dashed lines show the 
SHM, VL220 and VL270 distributions respectively. The left edge of the grey region corresponds to 
= 500 GeV. 

tion is greater for the VL220 distribution compared to the VL270 distribution, particularly 
in the high velocity region. However, the tail of the distribution is still much less populated 
in VL220 ^-iid VL270 compared to the SHM, which explains why the exclusion limits for 
VL220 and VL270 typically occur at larger cross sections than for the SHM. The increased 
population in the VL270 tail, relative to the VL220 tail, and the small change in particle 
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numbers over the region of Vmin of interest leads to stronger exclusion limits at high WIMP 
masses compared to VL22o- This shows that even a small change in the dispersion can lead 
to relatively large changes in allowed parameter space. 



VL270, = 80 GeV VL270, = 200 GeV 




S [keV] 6 [keV] 

Figure 6: Exclusion limits for VL27o. These results are similar to those presented in Figure 2. We 
see that there are still allowed regions of the DAMA parameter space at the 90% confidence level 
for WIMP masses lower than 200 GeV and S ^ 130 keV, however these typically occur at larger 
cross sections than in Figure 2. 
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Figure 7: The allowed parameter space for fixed S = 100 keV and S = 135 keV while varying M^. 
Comparing with Figure 4, we see that these results are similar to those obtained for the SHM, the 
difference being that the exclusion curves are typically at larger cross sections. While all of the 
DAMA parameter space is excluded at 5 = 100 keV (the right panel), there is clearly some allowed 
parameter space at ^ = 135 keV (the left panel) for WIMP masses less than 200 GeV. 
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5.3 The Dark Disc and iDM 

In Figure 8 we plot tlie exclusion limits as a function of S for the SHM and SHM + Dark 
Disc. As a fiducial choice we have taken R = £Smc_ _ 

PSHM 



SHM, = 90 GeV Dark Disc, = 90 GeV 




20 40 60 80 100 120 140 20 40 60 80 100 120 140 



S [keV] S [keV] 

Figure 8: Exclusion limits for the SHM (left panel) and the SHM + Dark Disc (right panel) for 
a WIMP mass of 90 GeV. The only major differences occur at relatively low 5 where one can see 
the region corresponding to DAM A channeled events is deformed when the Dark Disc is included. 
However, this region of parameter space is completely excluded by the other experiments. 




,5 [keV] 



Figure 9: The minimum velocity for scattering at CRESST-II (Dashed) and DAMA (SoUd). The 
velocities are calculated using Eq.(2.1) and taking the experimental energy threshold. The DAMA 
line corresponds to channeled scattering and the gray line is roughly the maximum velocity at which 
the Dark Disk particles would be found. For (5 > 15 keV very few particles in the dark disc will 
scatter. 



One can immediately see that the dark disc only influences the limits at very low 5. 
Increasing the density of the dark disc will change the limits at low 5 but not the value 
of 5 at which the disc is detectable. The fact that the dark disc is not detectable at high 
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6 is due to the increased minimum velocity required for scattering. As the disc lags the 
motion of the Sun by 40 km/s, and has velocity dispersion of 39 km/s in the tangential 
direction, 99% of disc particles will only have a relative velocity up to vdd ~ 135 km/s. 
If the minimum velocity for scattering is above this the dark disc will have little effect on 
event rates. This effect is demonstrated in Figure 9. 

One can see that for heavier WIMP masses the effect of the Dark Disc would increase, 
but would still only be found for 6 < 20 keV. These results would suggest that possible sub- 
halo components which are rotating roughly in the same way as the visible galactic plane 
would have little or no effect on inelastic scattering event rates, whereas the effect on elastic 
scattering rates should be pronounced. This effect could have interesting consequences for 
clumpy dark matter. 



5.4 Iodine quenching factor 

In [32] the iodine quenching factor is given as qj = 0.09 it 0.01. As changing the value 
of the quenching factor will stretch the recoil spectrum as a function of recoil energy it 
is interesting to see what effect varying the quenching factor has on the DAMA/LIBRA 
preferred region of parameter space, and whether this strongly influences the agreement 
with other experiments. The effect of this variation is illustrated in Figure 10. 
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Figure 10: The exclusion limits with an iodine quenching factor of qi = 0.08 (left panel) and 
qi = 0.1 (right panel), for a WIMP mass of 150 GeV and using the SHM. One can see that 
the preferred region of parameter space for the DAMA results can move by A ((5) ~ 30 keV and 
A(ct„) ^ 5 X 10~^''cm^, leading to agreement with or exclusion by the other experiments. 



As one can see, the overall effect can be significant, and the DAMA/LIBRA preferred 
region can shift by a relatively large amount. For the lower quenching factor there is a 
region of agreement between all experiments at the 90% confidence level, whereas for the 
higher quenching factor there is disagreement between DAMA/LIBRA and the CRESST-II 
and CDMS II experiments at almost a 99.5% level for the SHM. 
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As this effect is so pronounced it would be interesting if the DAMA collaboration 
are able to put tighter bounds on the quenching factors in their experiment. For elastic 
scattering the uncertainty is not so significant as it is the channeled event regions where qj = 
1 that show best agreement with the other experiments. However, as inelastic scattering 
finds greatest agreement for the region of parameter space corresponding to quenched 
events, the value of the quenching factor is of much greater importance. 

5.5 Circular velocity 

In [14] the Sun's circular velocity about the center of the Milky Way is given as Vdrc = 
220 lb 20km/s. In the iDM scenario the majority of observed events are expected to be 
from scattering off WIMPs at the high end of the tangential velocity distribution, so it is 
expected that varying the circular velocity will have an effect on the exclusion limits set by 
detectors. In Figure 11 we plot the exclusion curves for a WIMP of mass = 125 GeV 
for Vcirc = 200, 240km/s. 



SHM, = 125 GeV, Vcirc = 200 km/s SHM, = 125 GeV, v^rc = 240 km/s 




S [keV] S [keV] 



Figure 11: The exclusion limits for a a local circular velocity of Vdrc = 200 (left panel) and 
Vcirc = 220 (right panel), for a WIMP mass of 125 GeV and using the SHM. One can see that the 
limits change a small amount and that the overall agreement with the DAMA results is slightly 
better for a lower circular velocity. However, the overall qualitative effect of changing the circular 
velocity is relatively small. 

As expected this variation does affect the exclusion limits, however there is little overall 
change in the qualitative features. Unlike with the Dark Disc this variation occurs over 
all ranges of 5. Due to the enhancement of the signal, through increasing the number 
of particles with v > Vmin, increasing the circular velocity of the Sun generally leads to 
marginally stronger upper limits on cross sections for iDM. 

5.6 Local galactic escape velocity 

The RAVE survey have measured the escape velocity to lie within the range 498 < v^sc < 
608 km/s at 90% confidence [15]. In this paper, we have taken the fiducial value Vesc = 550 
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km/s, however the experiments are only sensitive to the tails of the velocity distribution, 
which depend on the escape velocity. In Figure 12 we illustrate what the effect of varying 
the escape velocity within the range allowed by the RAVE survey has on the experimental 
limits. 



SHM, (5 = 120 keV, v^,^ = 500 km/s 



SHM, 6= 120 keV, 
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Figure 12: Exclusion limits for the SHM for Ugsc = 500 km/s (left panel) and v^sc = 600 km/s 
(right panel). There is better agreement between experiments for lower escape velocity and the 
limits from CDMS II show the most sensitivity to the escape velocity. 



The limit which shows the largest change is from CDMS II, although all experimental 
limits show some variation, particularly at lower masses. This effect can be understood 
with reference to Figure 5. At lower masses, the minimum velocity is closer to the escape 
velocity, therefore for lower escape velocities, the range of integration is smaller, hence the 
limits are weaker. At higher masses, the relative change in varying the escape velocity 
is lower because the range of integration is larger, so the difference for all experiments is 
smaller. 



6. Conclusions 



We have shown, in the context of iDM, that the region of agreement between the DAMA 
data and results from other experiments is sensitive to the uncertainties present in the 
galactic WIMP velocity distribution. In particular we have found that the other direct de- 
tection experiments, including the most recent ZEPLIN-II and III, CRESST-II, XENONIO, 
KIMS, and CDMS-II data sets, do not exclude the region of parameter space preferred by 
the DAMA results up to WIMP masses ~ 1 TeV when the VL220 velocity distribution 
is used (see Figures 3 and 4), while for the VL270 velocity distribution, WIMP masses 
< 200 GeV are allowed (see Figures 6 and 7). Furthermore, we have also argued that 
the region of agreement between experiments is very sensitive to the quenching factor used 
to interpret the DAMA data (see Figure 10), and also to the local galactic escape velocity 
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(see Figure 12), both of which presently have ~10% uncertainties, and somewhat to the 
Sun's circular velocity (see Figure 11). 

Independent of experimental set-ups we have also shown that, in the iDM scenario, 
detectors would be insensitive to WIMPs traveling with a small velocity relative to the 
sun due to the increased minimum velocity for scattering. This would lead to dark matter 
clumps or streams being undetectable if rotating in the galactic plane, as in the Dark Disc 
halo model, or streaming with a small relative velocity (see Figure 8). 

New data from heavy element experiments such as CRESST-II and XENONIOO could 
have the potential to completely exclude the DAMA preferred region for iDM, especially if 
taken during the modulation maximum. In particular, if the planned EURECA experiment 
were to use a tungsten detector it would provide significant limits on both the elastic and 
iDM scattering cross-sections. Furthermore because of the low energy cut-off in the recoil 
energy spectrum for iDM, it is important that the experiments increase their sensitivity at 
higher recoil energies, where the iDM recoil spectrum peaks. 

Without a better understanding of the details of the velocity distribution in our galaxy, 
including the effects of the baryons, precise statements about the consistency of various 
direct detection experiments are not possible. 

Note added: Subsequent to the submission of this work two recent DM simulations were 
brought to our attention: 

• Recently the results of a second simulation. Via Lactea II [52], were made public. 
These results showed further detail in the substructure of the DM halo, with the 
dumpiness extending over six orders of magnitude, down to the resolution limit. 

• Another Milky Way-sized DM halo simulation, the Aquarius Project, also recently 
released their results [53], and the phase-space structure of these simulated halos 
was the subject of further study [54]. The velocity distributions found varied from 
simulation to simulation and showed interesting bumpy features. Further to this, 
in [54] the influence of these velocity distributions on direct detection signals for 
elastically scattering dark matter was investigated and it was found that the recoil 
spectra could deviate by up to 10% from that expected from the best-fit multivariate 
Gaussian. 

These simulations provide further strong support for deviations from the SHM, but as 
data from these simulations are not (to the authors knowledge) publicly available, we have 
not been able to update our study to include them, although we hope to return to this 
subject in a future publication. 
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